Código
# instalar/cargar paquetes
sketchy::load_packages(
c("ggplot2",
"viridis",
"caret"
)
)Loading required package: caretLoading required package: latticeSobreajuste y entrenamiento de modelos
Objetivo del manual
Entender el concepto de sobreajuste en modelos de aprendizaje estadístico.
Aprender a detectar sobreajuste en modelos de aprendizaje estadístico.
Aplicar metodos estadisticos para evitar el sobreajuste en modelos de aprendizaje estadístico.
Paquetes a utilizar en este manual:
Funciones personalizadas a utilizar en este manual:
Código
rmse_lm <- function(model, data, response = "y"){
pred <- predict(model, newdata = data)
sqrt(mean((data[, response] - pred)^2))
}
r2_lm <- function(model, data, response = "y"){
pred <- predict(model, newdata = data)
caret::R2(pred, data[, response])
}1 Entrenamiento de modelos
El objetivo de los modelos de aprendizaje estadístico es el de obtener patrones de los datos de entrenamiento para predecir o inferir correctamente los patrones en la población original de donde provienen esos datos de entrenamiento. Es decir, la clave esta en obtener patrones generales que sean extrapolables a nuevos datos. La idea principal del entrenamiento es ajustar el modelo a los datos de entrenamiento para aprender patrones que se puedan generalizar a datos nuevos. Sin embargo, parte de este proceso implica estrategias para evitar tanto el sobreajuste como el subajuste.
1.1 Sobreajuste
El sobreajuste ocurre cuando un modelo se ajusta demasiado bien a los datos de entrenamiento, capturando tanto los patrones verdaderos como el ruido o variaciones aleatorias de los datos. Como resultado, el modelo funciona bien en el conjunto de entrenamiento, pero tiene un rendimiento deficiente en nuevos datos (pobre capacidad de generalización). El sobreajuste se refiere a cuando modelo está tan ajustado a los datos de entrenamiento que afecta su capacidad de generalización. El sobreajuste se produce cuando un sistema de aprendizaje automático se entrena demasiado o con datos (levemente) sesgados, que hace que el algoritmo aprenda patrones que no son generales. Aprende características especificas pero no los patrones generales, el concepto.
Una forma de evaluar la capacidad de generalización de un modelo es mediante la división de los datos en dos conjuntos: entrenamiento y prueba. El modelo se ajusta a los datos de entrenamiento y se evalúa en los datos de prueba. El sobreajuste se puede detectar cuando el error en los datos de prueba es mucho mayor que el error en los datos de entrenamiento.
Los modelos más complejos tienden a sobreajustar más que lo modelos más simples. Además, ante un mismo modelo, a menor cantidad de datos es más posible que ese modelo se sobreajuste. Existen varios métodos para evaluar cuándo un modelo está sobreajustando. En la simulacion que se muestra a continuación, se ajusta un modelo de regresión lineal con diferentes cantidades de predictores (p). Se calcula el error cuadrático medio en los datos de entrenamiento y en los datos de prueba.
Código
repeticiones <- 100 # Número de repeticiones
n <- 100 # Número de observaciones
p <- 20 # Número de predictores
expr <- expression({
# Generar datos sintéticos
# Crear variables predictoras aleatorias
datos <- as.data.frame(matrix(rnorm(n * p), n, p))
colnames(datos) <- paste0("x", 1:p)
# Crear variable de respuesta con una combinación de algunas variables
datos$y <-
3 * datos$x1 - 2 * datos$x2 + 1 * datos$x3 + rnorm(n, 0, 2)
# Dividir en conjunto de entrenamiento y prueba
entren_indice <- createDataPartition(datos$y, p = 0.9, list = FALSE)
datos_entren <- datos[entren_indice, ]
datos_prueba <- datos[-entren_indice, ]
resultados_lista <- lapply(1:(ncol(datos) - 1), function(z) {
# Ajustar modelo de regresión lineal
modelo <-
lm(y ~ ., data = datos_entren[, c("y", paste0("x", 1:z))])
# Calcular raiz del error cuadrático medio
rmse_entren <- rmse_lm(modelo, datos_entren, response = "y")
rmse_prueba <- rmse_lm(modelo, datos_prueba, response = "y")
r2_entren <- r2_lm(modelo, datos_entren)
r2_prueba <- r2_lm(modelo, datos_prueba)
resultados <-
data.frame(
r2 = c(r2_prueba, r2_entren)
)
return(resultados)
})
resultados_df <- do.call(rbind, resultados_lista)
})
# repetir x veces
reps <- replicate(repeticiones, eval(expr), simplify = TRUE)
# promediar resultados
r2 <- apply(as.data.frame(reps), 1, mean)
resultados_promedio <-
data.frame(
n_predictores = rep(1:p, each = 2),
Tipo = c("Prueba", "Entrenamiendo"),
r2 = r2
)
ggplot(resultados_promedio, aes(x = n_predictores, y = r2, color = Tipo)) +
geom_line(lwd = 2) +
scale_x_continuous(breaks = seq(0, p, 1)) +
scale_color_viridis_d(end = 0.9) +
labs(x = "Número de predictores",
y = bquote('Coeficiente de\ndeterminación' ~ R ^ 2)) +
theme(
legend.background = element_rect(fill = "#fff3cd"),
legend.position = c(0.5, 0.2),
panel.background = element_rect(fill = "#fff3cd"),
plot.background = element_rect(fill = "#fff3cd", colour = NA)
) 
Podemos ver como en ambos casos el coeficiente de determinación (R2) aumenta en los primeros 3 predictores. Esto es de esperar ya que estos son los predictores asociados a la respuesta. Sin embargo, luego de este punto el R2 aumenta para los datos de entrenamiento, pero no para los datos de prueba. Esto es un claro indicio de sobreajuste.
En la siguiente simulación podemos ver mas claramente como el R2 calculado sobre datos de entrenamiento aumenta con la cantidad de predictores, a pesar de no haber un solo predictor asociado a la variable respuesta:
Código
repeticiones <- 100 # Número de repeticiones
n <- 100 # Número de observaciones
p <- 20 # Número de predictores
expr <- expression({
# Generar datos sintéticos
# Crear variables predictoras aleatorias
datos <- as.data.frame(matrix(rnorm(n * p), n, p))
colnames(datos) <- paste0("x", 1:p)
# Crear variable de respuesta con una combinación de algunas variables
datos$y <- rnorm(n, 0, 2)
# Dividir en conjunto de entrenamiento y prueba
entren_indice <- createDataPartition(datos$y, p = 0.9, list = FALSE)
datos_entren <- datos[entren_indice, ]
datos_prueba <- datos[-entren_indice, ]
resultados_lista <- lapply(1:(ncol(datos) - 1), function(z) {
# Ajustar modelo de regresión lineal
modelo <-
lm(y ~ ., data = datos_entren[, c("y", paste0("x", 1:z))])
# Calcular raiz del error cuadrático medio
rmse_entren <- rmse_lm(modelo, datos_entren, response = "y")
rmse_prueba <- rmse_lm(modelo, datos_prueba, response = "y")
r2_entren <- r2_lm(modelo, datos_entren)
r2_prueba <- r2_lm(modelo, datos_prueba)
resultados <-
data.frame(
r2 = c(r2_prueba, r2_entren)
)
return(resultados)
})
resultados_df <- do.call(rbind, resultados_lista)
})
# repetir x veces
reps <- replicate(repeticiones, eval(expr), simplify = TRUE)
# promediar resultados
r2 <- apply(as.data.frame(reps), 1, mean)
resultados_promedio_aleatorio <-
data.frame(
n_predictores = rep(1:p, each = 2),
Tipo = c("Prueba", "Entrenamiendo"),
r2 = r2
)
ggplot(resultados_promedio_aleatorio, aes(x = n_predictores, y = r2, color = Tipo)) +
geom_line(lwd = 2) +
scale_x_continuous(breaks = seq(0, p, 1)) +
scale_color_viridis_d(end = 0.9) +
labs(x = "Número de predictores",
y = bquote('Coeficiente de\ndeterminación' ~ R ^ 2)) +
theme(
legend.background = element_rect(fill = "#fff3cd"),
legend.position = c(0.7, 0.2),
panel.background = element_rect(fill = "#fff3cd"),
plot.background = element_rect(fill = "#fff3cd", colour = NA)
) 
Sin embargo, el R2 calculado sobre datos de prueba se mantiene constante.
2 AIC
El criterio de información de Akaike (AIC) es una medida de la calidad relativa de un modelo estadístico para un conjunto dado de datos. Dado un conjunto de modelos, el modelo preferido es el que tiene el menor AIC. Dado un conjunto de modelos, el AIC proporciona una medida relativa de la calidad de cada modelo, en relación con los otros modelos. Un AIC más bajo indica un mejor modelo.
Podemos ver como el AIC aumenta con la cantidad de predictores en el modelo. En el siguiente gráfico se muestra el AIC para un modelo de regresión lineal con diferentes cantidades de predictores:
Código
repeticiones <- 100 # Número de repeticiones
n <- 100 # Número de observaciones
p <- 20 # Número de predictores
expr <- expression({
# Generar datos sintéticos
# Crear variables predictoras aleatorias
datos <- as.data.frame(matrix(rnorm(n * p), n, p))
colnames(datos) <- paste0("x", 1:p)
# Crear variable de respuesta con una combinación de algunas variables
datos$y <-
3 * datos$x1 - 2 * datos$x2 + 1 * datos$x3 + rnorm(n, 0, 2)
# Dividir en conjunto de entrenamiento y prueba
entren_indice <- createDataPartition(datos$y, p = 0.9, list = FALSE)
datos_entren <- datos[entren_indice, ]
datos_prueba <- datos[-entren_indice, ]
resultados_lista <- lapply(1:(ncol(datos) - 1), function(z) {
# Ajustar modelo de regresión lineal
modelo <-
lm(y ~ ., data = datos[, c("y", paste0("x", 1:z))])
# modelos nulos
modelo_nulo <- lm(y ~ 1, data = datos[, c("y", paste0("x", 1:z))])
# Calcular AIC
aics <- AIC(modelo, modelo_nulo)
aics$delta.aic <- aics$AIC - min(aics$AIC)
resultados <- aics["modelo_nulo", "delta.aic"]
return(resultados)
})
resultados_df <- do.call(rbind, resultados_lista)
})
# repetir x veces
reps <- replicate(repeticiones, eval(expr), simplify = TRUE)
# promediar resultados
delta_aic <- apply(as.data.frame(reps), 1, mean)
resultados_promedio_aic <-
data.frame(
n_predictores = rep(1:p),
delta_aic = delta_aic
)
ggplot(resultados_promedio_aic, aes(x = n_predictores, y = delta_aic)) +
geom_line(lwd = 2, color = viridis(1)) +
scale_x_continuous(breaks = seq(0, p, 1)) +
labs(x = "Número de predictores", y = "Delta AIC") +
theme(legend.position = c(0.5, 0.7)) +
scale_y_reverse()
Datos donde ningun predictor está asociado:
Código
repeticiones <- 100 # Número de repeticiones
n <- 100 # Número de observaciones
p <- 20 # Número de predictores
expr <- expression({
# Generar datos sintéticos
# Crear variables predictoras aleatorias
datos <- as.data.frame(matrix(rnorm(n * p), n, p))
colnames(datos) <- paste0("x", 1:p)
# Crear variable de respuesta no asociada a los predictores
datos$y <-rnorm(n, 0, 2)
# Dividir en conjunto de entrenamiento y prueba
entren_indice <- createDataPartition(datos$y, p = 0.9, list = FALSE)
datos_entren <- datos[entren_indice, ]
datos_prueba <- datos[-entren_indice, ]
resultados_lista <- lapply(1:(ncol(datos) - 1), function(z) {
# Ajustar modelo de regresión lineal
modelo <-
lm(y ~ ., data = datos[, c("y", paste0("x", 1:z))])
# modelos nulos
modelo_nulo <- lm(y ~ 1, data = datos[, c("y", paste0("x", 1:z))])
# Calcular AIC
aics <- AIC(modelo, modelo_nulo)
aics$delta.aic <- aics$AIC - min(aics$AIC)
resultados <- aics["modelo_nulo", "delta.aic"]
return(resultados)
})
resultados_df <- do.call(rbind, resultados_lista)
})
# repetir x veces
reps <- replicate(repeticiones, eval(expr), simplify = TRUE)
# promediar resultados
delta_aic <- apply(as.data.frame(reps), 1, mean)
resultados_promedio_aic <-
data.frame(
n_predictores = rep(1:p),
delta_aic = delta_aic
)
ggplot(resultados_promedio_aic, aes(x = n_predictores, y = delta_aic)) +
geom_line(lwd = 2, color = viridis(1)) +
scale_x_continuous(breaks = seq(0, p, 1)) +
labs(x = "Número de predictores", y = "Delta AIC") +
theme(legend.position = c(0.5, 0.7)) +
scale_y_reverse(limits = c(1, -2))
2.1 Ejercicio 1
1.1
Información de la sesión
R version 4.3.2 (2023-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 22.04.2 LTS
Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
time zone: America/Costa_Rica
tzcode source: system (glibc)
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] caret_6.0-94 lattice_0.20-45 nnet_7.3-17 viridis_0.6.5
[5] viridisLite_0.4.2 ggplot2_3.5.1 knitr_1.48
loaded via a namespace (and not attached):
[1] gtable_0.3.5 xfun_0.47 htmlwidgets_1.6.4
[4] recipes_1.0.9 remotes_2.5.0 vctrs_0.6.5
[7] tools_4.3.2 generics_0.1.3 stats4_4.3.2
[10] parallel_4.3.2 tibble_3.2.1 fansi_1.0.6
[13] ModelMetrics_1.2.2.2 pkgconfig_2.0.3 Matrix_1.6-5
[16] data.table_1.14.10 lifecycle_1.0.4 farver_2.1.2
[19] compiler_4.3.2 stringr_1.5.1 munsell_0.5.1
[22] codetools_0.2-18 sketchy_1.0.3 htmltools_0.5.8.1
[25] class_7.3-20 yaml_2.3.10 prodlim_2023.08.28
[28] pillar_1.9.0 crayon_1.5.3 MASS_7.3-55
[31] gower_1.0.1 iterators_1.0.14 rpart_4.1.16
[34] foreach_1.5.2 nlme_3.1-155 parallelly_1.38.0
[37] lava_1.7.3 tidyselect_1.2.1 packrat_0.9.2
[40] digest_0.6.37 stringi_1.8.4 future_1.34.0
[43] reshape2_1.4.4 purrr_1.0.2 dplyr_1.1.4
[46] listenv_0.9.1 labeling_0.4.3 splines_4.3.2
[49] fastmap_1.2.0 grid_4.3.2 colorspace_2.1-1
[52] cli_3.6.3 magrittr_2.0.3 survival_3.2-13
[55] utf8_1.2.4 future.apply_1.11.2 withr_3.0.1
[58] scales_1.3.0 timechange_0.2.0 lubridate_1.9.3
[61] rmarkdown_2.28 globals_0.16.3 timeDate_4032.109
[64] gridExtra_2.3 evaluate_1.0.0 hardhat_1.3.0
[67] rlang_1.1.4 Rcpp_1.0.13 glue_1.8.0
[70] xaringanExtra_0.8.0 pROC_1.18.5 ipred_0.9-14
[73] rstudioapi_0.16.0 jsonlite_1.8.9 R6_2.5.1
[76] plyr_1.8.9